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LONG-TERM  GOALS 

The  over-all  goal  of  this  research  is  the  development  of  an  accurate  and  reliable  propagation  model 
applicable  to  environments  which  exhibit  strong  range  dependence  in  all  three  spatial  dimensions. 

OBJECTIVES 

The  objective  of  this  work  is  to  gain  an  understanding  the  physics  of  propagation  in  continental  shelf 
areas,  specifically  horizontal  refraction  and  mode  coupling  induced  by  three-dimensional  (3D) 
inhomogeneities  in  the  waveguide.  A  coupled-mode  approach  has  been  applied  for  this  purpose.  The 
coupled-mode  approach  is  attractive  for  solving  problems  involving  3D  propagation  for  several  reasons. 
First,  this  technique  provides  intuitive  results  for  understanding  the  features  responsible  for  observed 
propagation  effects  in  range-dependent  environments.  For  example,  upslope  propagation  is 
characterized  by  acoustic  energy  radiated  into  the  bottom  at  discrete  depths  associated  with  mode 
cut-off.  The  modal  decomposition  of  the  acoustic  field  has  also  been  used  to  describe  horizontal 
refraction  in  a  wedge-shaped  ocean,  for  which  the  single-mode  interference  pattern  associated  with  rays 
launched  up  and  across  the  shelf  has  been  well  documented  [Weinberg  and  Burridge  (1974)]. 
Furthermore,  coupled-mode  solutions  can  be  highly  accurate  and  have  been  used  for  benchmarking 
solutions  to  range-dependent  problems  [Jensen  and  Ferla  (1990)].  In  order  to  appreciate  the  limitations 
of  existing  3D  models,  it  is  necessary  to  have  methods  which  can  provide  reference  solutions  for 
comparison. 

APPROACH 

A  3D  acoustic  propagation  model  based  on  the  stepwise  coupled-mode  approach  [Evans  (1983)] 
implemented  as  a  single-scatter  solution  has  been  developed.  This  technique  is  based  on  a  hybrid 
modeling  approach  for  which  normal  modes  supply  the  vertical  dependence  and  a  parabolic  equation 
(PE)  solution  provides  the  horizontal  dependence.  The  inhomogeneous  Helmholtz  equation  for  pressure 
P(r.  0,z)  at  range  r,  azimuth  0,  and  depth  z  from  a  point  continuous  wave  source  of  amplitude  5(a)), 
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located  at  range  r  =  0  and  depth  z  —  zo,  is  given  by 


p(r,0,z)V- 
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p{r,e,z) 


yp(r,e,z) 


+  k2(r,e,z)P(r,G,z)  =  -4^S(©)^^5(z-zo),  (1) 


where  V  =  r  +  zjp  is  the  gradient  operator  in  cylindrical  coordinates,  k(r.  6,z)  —  0)/c(r.  0,z.) 

is  the  acoustic  wavenumber,  co  =  2 nf,  f  is  the  acoustic  frequency,  c(r,  O.z.)  is  sound  speed,  and 
p(r,  O.z.)  is  density.  In  this  work,  the  range-dependent  environment  is  approximated  by  a  number  of 
range-independent  regions  defined  in  the  range-bearing  plane.  In  the  discrete  environment,  the 
continuous  variables  of  Eq.  (1)  are  discretized  in  range  and  bearing  such  that  the  range-independent 
region  at  the  ith  step  in  range  and  /h  step  in  azimuth  is  characterized  by  sound  speed  Cij(z ),  density 
Pij(z),  and  wavenumber  kjj(z).  The  solution  for  pressure  in  the  zth  step  in  range  is  found  by  a 
separation  of  variables 


M 

Pi(r,  0,z)  =  £  Amj (r ,  Q)tymUj{z)  for  r,-_ i  <  r  <  ru 

m=\ 


(2) 


where  Anli(r,  0)  are  the  modal  amplitudes  in  the  zth  step  in  range,  and  <j>mij(z)  are  the  modal 
eigenfunctions  in  the  zth  step  in  range  and  /h  step  in  azimuth.  As  denoted  by  Eq.  (2),  the  modal 
amplitudes  are  continuous  functions  of  range  and  azimuth  whereas  modal  eigenfunctions  are 
continuous  functions  in  depth  only. 


The  modal  eigenfunctions  (f>mj ;  (z)  are  defined  locally  in  each  range-independent  region  such  that  they 
satisfy 


piAz)Tz 


r  1  d(j)mij(zY 
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(3) 


where  kmjj  is  the  horizontal  wavenumber  of  the  mth  mode.  Within  each  region,  the  environment  consists 
of  an  arbitrary  number  of  fluid  layers  bounded  above  by  a  pressure  release  surface  and  below  by  a  fluid 
half-space.  In  this  work,  the  Pekeris  branch  cut  is  chosen  such  that  the  total  pressure  is  calculated  from 
a  Pekeris  branch  line  integral,  plus  a  finite  sum  of  trapped  modes,  plus  an  infinite  sum  of  leaky  modes. 

A  small  gradient  is  introduced  in  the  lower  half-space  which  effectively  removes  the  branch  point  and 
associated  branch  cut  from  the  problem  [Westwood  and  Koch  (1999)]. 


In  the  solution  for  the  modal  amplitudes,  environmental  range  dependence  takes  place  discretely  at  the 
interfaces  between  regions.  Within  the  zth  step  in  range,  the  environment  varies  as  a  function  of  azimuth 
only,  and  the  modal  amplitudes  satisfy  the  2D  Helmholtz  equation 
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A 0)  =  0  for  n-i  <  r  <  r,. 


(4) 


This  equation  must  be  solved  for  each  mode  with  the  horizontal  refraction  determined  by  the  modal 
phase  speed  cf^j  —  ft)/ Re { knv  j }  and  modal  attenuation  amjj  —  lm {k„Uj}.  The  solution  is  obtained  from 
a  PE  model  [Collins  (1994)]  by  factoring  Eq.  (4)  and  keeping  the  outgoing  factor  to  obtain 


dA'njr  ^  =  ik0 \j  1  +k()2(J^-^+k2  -/c/) Amj  ( r,  0 )  for  r,-_  i  <  r  <  rt ,  (5 ) 

where  ko  —  —  is  the  reference  wavenumber,  and  cq  is  the  reference  sound  speed.  Equation  (5)  differs 

H) 

from  the  original  PE  model  which  is  defined  in  the  range-depth  plane  with  pressure  release  boundary 
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conditions  at  the  top  and  bottom  of  the  waveguide.  The  starting  function  is  a  delta  function  located  at 
the  origin  with  amplitude  determined  by  the  source  depth  and  the  modal  eigenfunctions  in  region  i  —  1, 


1W1  (0,0)  =4?rS(a>) 

PijUoJ 


(6) 


In  the  stepwise  coupled-mode  model,  mode  coupling  is  calculated  in  the  outgoing  direction  at  the 
interfaces  between  the  range-independent  regions.  Azimuthal  mode  coupling  is  neglected,  such  that  the 
mode-coupling  occurs  only  through  the  boundary  conditions  at  interfaces  defined  in  the  tangential 
direction.  The  two-way  coupled-mode  solution  is  obtained  by  requiring  continuity  of  pressure  and 
particle  velocity  at  r  =  rj  located  at  the  boundary  between  the  (/,  j)  and  (/  +l,j)  range-independent 
regions.  In  this  work,  the  backscattered  field  is  neglected  to  obtain  the  single-scatter  solution 
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WORK  COMPLETED 


The  main  accomplishments  of  2013  include:  (1)  a  comparison  of  the  3D  coupled-mode  solution  to  a 
finite  element  method  (FEM)  solution,  (2)  ranging  estimates  from  horizontal  multipath  in  a  continental 
shelf  environment,  (3)  modeling  of  3D  propagation  effects  observed  in  the  Church  Stroke  III  data  set, 
and  (4)  development  of  a  robust  normal  mode  code  for  environments  with  elastic  layers. 


RESULTS 


Comparison  of  the  3D  coupled-mode  solution  to  a  FEM  solution 

As  described  above,  in  formulating  the  3D  coupled-mode  solution,  it  was  assumed  that  mode-coupling 
in  the  azimuthal  direction  is  small  and  can  be  neglected.  The  effects  of  this  assumption  are  examined  by 
comparing  the  solution  from  the  3D  coupled-mode  model  to  a  solution  obtained  using  a  FEM  model 
applied  with  the  longitudinally  invariant  (LI)  technique.  This  technique  is  appropriate  for  environments 
which  have  lateral  variations  that  occur  in  one  dimension  only.  The  solution  is  obtained  by  applying  a 
cosine  transform  to  eliminate  the  range-independent  dimension  [Isakson  and  Chotiros  (201 1); 
Goldsberry  and  Isakson  (2013)]. 

The  propagation  environment  consists  of  an  isovelocity  water  column  (c  =  1500  m/s,  p  =  1.0  g/cm3) 
over  an  acoustic  half-space  having  the  properties  of  sand  (c  =  1767.3  m/s,  p  =  1.845  g/cm3,  a  = 

0.75225  dB/A).  The  only  range-dependent  feature  is  a  cosine-shaped  hill  that  extends  infinitely  in  the 
y-direction,  such  that  the  water  depth  d(x)  is  defined  by 


d{x) 


47.5-12.5cos(§^)  if  2750  <x<  3250 
60  if  x  <  2750  or  ^  >  3250. 


(8) 


The  solution  was  calculated  for  a  source  located  at  a  depth  of  20  m  and  a  frequency  of  100  Hz.  The 
infinite  cosine  hill  induces  strong  horizontal  refraction  and  mode-coupling  effects.  This  is  demonstrated 
by  comparisons  to  Nx2D  coupled-mode  and  3D  adiabatic-mode  solutions  in  Fig.  1(a)  and  (b), 
respectively. 
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(a)  Comparison  to  Nx2D  Coupled-Mode  Solution 


(b)  Comparison  to  3D  Adiabatic-Mode  Solution 


(c)  Comparison  to  LI  FEM  Solution 


Figure  1:  Comparison  of  the  3D  coupled-mode  solution  at  a  depth  of  20  m  along  x  -  5  km  to 
solutions  calculated  by  (a)  Nx2D  adiabatic-mode  technique  and  by  (b)  3D  adiabatic-mode  technique 
to  demonstrate  the  relative  importance  of  these  effects.  Comparison  of  the  3D  coupled-mode  solution 
to  a  solution  calculted  using  the  FEM  with  LI  technique  to  demonstrate  solution  accuracy. 


The  3D  coupled-mode  solution  for  transmission  loss  (TL)  over  the  infinite  cosine  hill  at  a  depth  of  20  m 
along  x  =  5  km  is  compared  to  the  solution  calculated  by  the  LI  technique  with  the  FEM  in  Fig.  1(c). 
With  the  exception  of  the  nulls  in  the  TL  curves,  the  solutions  are  within  0.5  dB.  This  error  is  on  the 
same  order  as  the  error  in  2D  solutions  calculated  for  propagation  over  the  axisymmetric  cosine  hill. 
Therefore,  the  error  in  the  3D  solutions  is  attributed  to  the  intrinsic  differences  in  the  approaches  used  to 
solve  the  2D  propagation  problem.  Hence,  the  assumption  of  negligible  azimuthal  mode  coupling  was 
validated  for  this  environment. 

Ranging  estimates  from  horizontal  multipath  in  a  continental  shelf  environment 

This  work  seeks  to  exploit  horizontal  multipath  effects  measured  in  beamformed  data  from  a  horizontal 
line  array  (HLA)  to  estimate  the  range  of  an  acoustic  source.  A  normal  mode  approach  is  applied,  and 
rays  are  traced  in  the  horizontal  plane  with  refraction  determined  by  the  modal  phase  speed.  The 
technique  is  applied  to  acoustic  data  recorded  on  22  February  2009  by  an  HLA  located  12  km  east  of 
the  southern  coast  of  Florida.  The  array  was  oriented  downslope,  providing  a  broadside  look  along  the 
250  m  isobath.  During  the  experiment,  the  RJV  Walton  Smith  towed  a  source  broadcasting  a  2-minute 
sequence  of  signals  [Heaney  and  Murray  (2009)].  The  source  was  towed  at  a  depth  of  100  m  as  the  R/V 
Walton  Smith  traveled  north-northwest,  parallel  to  the  Florida  coastline  until  it  was  approximately  63 
km  north  of  the  array.  The  ship  then  turned  around  and  towed  the  source  at  a  depth  of  20  m  as  it 
traveled  back  south- southeast  along  the  same  path. 
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Figure  2:  (a)  BTR  with  bearing  angles  of  the  direct  and  refracted  paths  highlighted  by  the  dashed 
lines.  The  horizontal  dashed  line  at  time  08:38  indicates  the  data  used  for  additional  analysis,  (b) 
Horizontal  ray  traces  for  modes  one  through  four  at  50  Hz  overlaid  on  bathymetry.  The  jc-  and  y-axes 
are  distance  relative  to  the  array  location,  (c)  Comparison  of  estimated  source  range  using  the 
multipath  ranging  technique  to  calculated  range  using  GPS  data.  The  shaded  region  represents  10% 

ranging  error  referenced  to  the  GPS  data. 


A  bearing  time  record  (BTR)  of  the  data  is  displayed  in  Fig.  2(a).  The  bearing  angles  of  the  direct  and 
refracted  paths  of  the  signal  broadcasted  by  the  towed  source  are  highlighted  by  the  dashed  lines.  To 
illustrate  the  multipath  ranging  technique,  the  horizontal  ray  trace  for  a  single  time  scan  is  shown  in 
Fig.  2(b).  The  time  selected  for  this  analysis  is  08:38  and  is  indicated  by  the  horizontal  dashed  line  in 
Fig.  2(a).  At  this  time,  the  estimated  bearing  angles  of  the  direct  and  refracted  paths  were  86.63°  and 
78.57°,  respectively.  The  ray  trace  was  performed  using  the  modal  phase  speed  for  modes  one  through 
four  at  50  Hz  to  determine  the  refraction  of  the  rays.  The  rays  traced  for  the  direct  and  refracted  paths 
are  shown  by  the  solid  and  dashed  lines,  respectively,  for  all  four  modes.  The  location  of  the  source  is 
estimated  from  the  intersection  of  the  rays,  marked  by  the  “x’s”  in  Fig.  2(b). 

The  procedure  described  above  was  repeated  for  60  estimates  of  the  direct  and  refracted  bearing  angles 
over  an  eight  hour  period  to  build  up  a  track  of  the  source’s  location  over  time.  For  this  application,  the 
modal  phase  speed  corresponding  to  mode  three  at  50  Hz  was  used.  The  ranging  estimate  is  presented 
in  Fig.  2(c).  The  largest  errors  in  the  estimate  occur  when  the  source  is  farthest  away  from  the  array  and 
lower  SNR  has  degraded  the  accuracy  of  the  estimated  bearing  angles.  For  this  result,  the  average  of  the 
absolute  error  in  the  range  estimate  is  3.7%,  and  the  error  does  not  exceed  10%  at  any  range. 
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(a)  Bathymetry 


(d)  Modeled  Data 
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(b)  Stratigraphy 


1000 

^N, 

Bathymetry 

1500 

Subbottom 

r  " 

2000 

2500 

jE 

£  3000 
ft 

/ 

°  3500 

/ 1 

4000 

Ooze 

4500 

Dolomite 

5000 

-40  -20  0  20  40 

Range  [km] 


10  20  30  40  50 

Time  [s] 


Time  [s] 


Figure  3:  (a)  Bathymetry  of  the  Catoche  Tongue  with  the  location  of  the  AC OD AC  array  (triangle) 
and  SUS  charges  (dots).  ( b )  Reconstructed  stratigraphy  of  the  Catoche  Tongue  along  a  northwest 
transect  passing  through  the  ACODAC  array  location  with  the  array  located  at  0  km  and  positive 
ranges  oriented  to  the  southeast,  (c)  Data  recorded  on  channel  11  of  the  ACODAC  array  from  the  13 
SUS  charges  deployed  at  the  locations  indicated  in  Fig.  3(a).  The  y-axis  is  distance  from  the  SUS 
charge  to  the  ACODAC  array,  (d)  Modeled  data  for  the  SUS  charge  located  57  km  away  from  the 
ACODAC  array  beamformed  with  a  HLA  orientated  with  90°  in  the  northwest  direction. 


Modeling  of  3D  propagation  effects  observed  in  the  Church  Stroke  III  data  set 

Three  dimensional  horizontally  refracted  arrivals  recorded  in  the  Catoche  Tongue  region  of  the  Gulf  of 
Mexico  as  part  of  the  Church  Stroke  III  Exercise  in  1979  were  modeled.  Acoustic  data  from  the 
Acoustic  Data  Capsule  (ACODAC)  vertical  line  array  (VLA),  located  in  the  center  of  the  tongue,  as 
shown  in  Fig.  3(a),  were  recently  recovered  as  part  of  the  Long-Range  Acoustic  Propagation  Project 
(LRAPP)  tape  recovery  effort.  The  ACODAC  array  was  moored  at  23.62°  N,  86.02°  W  in 
approximately  3470  m  of  water.  During  the  exercise,  361  Signal  Underwater  Sound  (SUS)  charges  were 
deployed.  A  subset  of  13  shots  located  nearest  the  array  were  examined  in  this  work,  and  their  locations 
are  indicated  in  Fig.  3(a).  Sagers  and  Ballard  (2013)  assessed  the  stratigraphic  structure  of  the  seabed 
was  from  seismograms  and  core  data  [Buffler  et  al.  (1984);  Worzel  (1973);  Shaub  (1987);  Brunson  et  al. 
(1980)].  Fig.  3(b)  shows  the  geoacoustic  representation  used  for  the  propagation  modeling. 

The  time  series  for  the  SUS  events  recorded  on  channel  1 1  (945  m  depth)  are  shown  in  Fig.  3(c).  The 
first  arrival  for  each  event  has  been  aligned  at  10  s.  The  direct  (non-horizontally  refracted)  arrivals  are 
received  within  the  first  10  seconds  after  the  first  arrival.  A  later  set  of  horizontally  refracted  arrivals  are 
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also  present  in  every  SUS  event  with  the  arrival  time  increasing  as  the  SUS  range  decreases.  The  later 
arrivals  also  appear  to  become  more  diffuse  as  the  SUS  range  decreases. 

The  acoustic  data  were  modeled  in  a  5  Hz  band  centered  at  25  Hz  using  a  hybrid  3D  adiabatic-mode 
approach  [Ballard  (2012)].  The  modeled  data  for  the  SUS  charge  located  57  km  from  the  ACODAC 
array  (indicated  by  the  arrow  in  Fig.  3(a))  are  shown  in  Fig.  3(d).  The  modeled  data  were  beamformed 
with  a  HLA  orientated  with  broadside  look  in  the  southwest  direction.  The  beamformed  data  clearly 
show  the  direct  path  arrival  at  angle  of  approximately  5°  with  a  strong  refracted  arrival  at  later  times 
between  -60°  and  -20°  and  a  weak  refracted  arrival  at  approximately  40°.  The  strong  refracted  arrival  is 
caused  by  interaction  with  the  southern  slope  of  the  Catoche  Tongue  while  the  weak  refracted  arrival 
results  from  interaction  with  the  northern  slope.  The  relative  strength  of  the  refracted  arrivals  is  caused 
by  the  comparative  steepness  of  the  southern  slope  and  by  its  proximity  to  the  SUS  charges. 

Development  of  normal  mode  code  for  environments  with  elastic  layers 


A  robust  normal  mode  code  for  calculating  the  eigenvalues  and  eigenfunctions  in  environments  with 
elastic  layers  using  a  FEM  approach  is  under  development.  In  this  solution,  the  environment  is  assumed 
to  be  made  up  of  a  number  of  layers  with  constant  density  p  and  shear  wave  speed  cs.  The 
compressional  wave  speed  cp  is  allowed  to  vary  within  a  layer.  Within  each  layer,  the  shear  y/  and 
compressional  (j)  potentials  satisfy 


V20 


v2w=-^ 

C2  dt2  ’  V  C2  dt2  ' 


(9) 


The  displacements  and  stresses  can  be  written  in  terms  of  the  potentials  as 


d<$>  d2y 
"(r'z)=a7  +  cW? 

azz(r,z)  =  (x+2^df+xl-d-M, 

dz  r  dr 

,  .  ,  Id,.  dw  „  du 

Orr{r,z)  =  X  -^-(ru)  +  —  +  2/^——, 

r  dr  dz  dr 


.  .  d(b  1  d  dw 

wf,Z)  =  -3 - , 

dz  r  dr  dr 

,  .  du  dw 

<Jrz{r,Z)=P  3- +  3-  , 

dz  dr 

\  1  1  d  ,  .  dw  u 

r-,z)  =  A  -^-(rw)  +  —  +2/I-. 

r  dr  dz  r 


where  A  and  p  are  the  Lame  constants. 


Solutions  to  range-dependent  problems  will  be  calculated  using  the  stepwise  coupled-mode  approach 
described  above  with  boundary  conditions  defined  by  the  momentum  conservation  equations  [Koch 
et  al.  (1983)] 


d2u 

'd? 

d2w 
'  ~dt2 


ld_ 
r  dr 


(r<J,T)  + 


ld_ 
r  dr 


(rOrz)  T 


dorz 

dz 

d2C 7ZZ 

dz2 


1 


IMPACT/APPLICATIONS 

The  impact  of  this  work  is  an  increased  understanding  of  acoustic  propagation  through  complicated 
coastal  environments  for  which  the  bathymetry,  seabed  properties,  and  oceanography  can  vary  in  three 
dimensions. 
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TRANSITIONS 


The  primary  transition  for  this  project  is  an  accurate  and  reliable  model  for  acoustic  propagation  in 
environments  with  strong  three-dimensional  range  dependence.  Because  coupled-mode  approaches  are 
computationally  intensive,  they  have  historically  been  used  to  benchmark  faster  techniques  which 
approximate  the  solution  to  the  wave  equation. 

RELATED  PROJECTS 

Geoacoustic  inversion  in  three-dimensional  environments 

The  goal  of  this  project  is  to  estimate  water  column  sound  speed  in  a  3D  volume  using  modal  travel 
time  measurements  from  multiple  source-receiver  pairs.  A  thorough  understanding  of  the  forward 
problem,  including  the  effects  of  horizontal  refraction  and  mode  coupling,  is  necessary  to  successfully 
estimate  environmental  parameters  in  regions  with  3D  inhomogeneities.  This  project  is  currently 
funded  by  ARL:UT’s  IR&D  program. 

Acoustic  propagation  modeling  for  diver  detection  sonar  systems 

The  purpose  of  this  work  is  to  characterize  waveforms  at  virtual  receiver  distances  on  the  order  of  a 
1000  meters  away  from  active  diver  detection  sonar  systems  installed  at  fixed  locations  within 
operational  sites  of  interest.  The  3D  coupled-mode  model  is  applied  for  this  purpose.  This  project  is 
funded  by  Space  and  Naval  Warfare  Systems  Center  Pacific. 
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